{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<center>\n",
    "<img src=\"../../img/ods_stickers.jpg\">\n",
    "## Открытый курс по машинному обучению. Сессия № 2\n",
    "# <center>Классификация групповых и одиночных целей\n",
    "### <center> Автор: Мустаев Айрат (@airat)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## <center> Индивидуальный проект по анализу данных </center>"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "from matplotlib import pyplot as plt\n",
    "import scipy.stats\n",
    "from scipy import signal\n",
    "from scipy.signal import butter, lfilter, freqz\n",
    "from numpy.fft import irfft,rfft, rfftfreq \n",
    "%matplotlib inline\n",
    "import scipy.stats as sts\n",
    "from pylab import *\n",
    "from scipy import fft\n",
    "from scipy import signal\n",
    "import copy\n",
    "from scipy.signal import hilbert, chirp"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# -*- coding: utf-8 -*-\n",
    "%pylab inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import copy\n",
    "import seaborn as sns;\n",
    "import itertools "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.tree import DecisionTreeClassifier # Ваш код здесь\n",
    "from sklearn.ensemble import RandomForestClassifier\n",
    "from sklearn.ensemble import AdaBoostClassifier\n",
    "from sklearn.ensemble import GradientBoostingClassifier\n",
    "from sklearn.ensemble import ExtraTreesClassifier\n",
    "from sklearn.cross_validation import cross_val_score\n",
    "from sklearn.grid_search import GridSearchCV"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.metrics import accuracy_score\n",
    "from sklearn.metrics import roc_auc_score\n",
    "from sklearn.metrics import confusion_matrix, f1_score"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.learning_curve import learning_curve\n",
    "from sklearn.learning_curve import validation_curve"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "В проекте рассмотрим решение задачи классификации двух классов:\n",
    " * одиночный человек;\n",
    " * группа людей.\n",
    " \n",
    "Исходными данными являются эксперементальные (сырые) данные, полученные с АЦП [(АЦП ЦАП Zet 230)](https://zetlab.com/shop/izmeritelnoe-oborudovanie/moduli-atsp-tsap/atsp-tsap-zet-230).\n",
    "\n",
    "Под сырыми данными, понимаются синхронизированные по времени реализация (сейсмического сигнала) с частотой дискретизации Fs = 500 Гц. \n",
    "\n",
    "### Краткое описание природы данных\n",
    "\n",
    "#### Сейсмическая волна\n",
    "Сейсмическая волна - это волны, переносящие энергию упругих (механических) колебаний в горных породах. Источником сейсмической волны может быть землетрясение, взрыв, вибрация или <b>удар (в нашем случае проход объекта классификации)</b>. \n",
    "<p>Существует следующая классификация сейсмический волн:</p>\n",
    " * Объёмные волны - Объёмные волны проходят через недра Земли. Путь волн преломляется различной плотностью и жёсткостью подземных пород.\n",
    " * P-волны (первичные волны) — продольные, или компрессионные волны. Обычно их скорость в два раза быстрее S-волн, проходить они могут через любые материалы.\n",
    " * P- и S-волны в мантии и ядре.\n",
    " * Поверхностные волны несколько похожи на волны воды, но в отличие от них они путешествуют по земной поверхности. Их обычная скорость значительно ниже скорости волн тела. Из-за своей низкой частоты, времени действия и большой амплитуды они являются самыми разрушительными изо всех типов сейсмических волн."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Загрузим и рассмотрим исходные данные (сигналы) и их метаданные:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "file_obj1 = open('file_signal/sig0002_20_1.txt', 'r')\n",
    "data_s = file_obj1.read().split(\"\\n\")\n",
    "\n",
    "Fs = 500.0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "meta_data = pd.read_csv('metaTable_erase.csv', sep = ';')\n",
    "meta_data.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "где:\n",
    " * nameFile - название файла в котором закодирована дата и время записи сигнала;\n",
    " * Count - количество людей для проводилась запись сигнала;\n",
    " * Steps - количество проходов в рамках одной записи сигнала."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "seism = {}  \n",
    "\n",
    "for i in range(len(data_s)):\n",
    "    sig_s = data_s[i].split('\\t')[1]\n",
    "    sig_s = list(map(lambda value: \\\n",
    "                float(value), \n",
    "                filter(lambda value2: value2 != '' and len(value2) > 0, \n",
    "                sig_s.split(';')\n",
    "                )\n",
    "    ))\n",
    "    seism[data_s[i].split('\\t')[0].upper()] = sig_s\n",
    "\n",
    "name_signal = seism.keys()\n",
    "name_signal = list(name_signal)\n",
    "print('Количество реализаций - {}'.format(len(name_signal)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Рассмотрим реализацию сигнала, например для реализации = S140314_161345:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plotSignal(nameEx):\n",
    "    seismic_signal = seism[nameEx]\n",
    "    time_x = np.linspace(0,len(seismic_signal)/Fs, len(seismic_signal))\n",
    "    plt.figure(figsize = (20, 6))\n",
    "    plt.title('Реализация сейсмического сигнала')\n",
    "    plt.plot(time_x,seismic_signal, )\n",
    "    plt.legend(['Сейсмический сигнал'])\n",
    "    plt.xlabel('Время, с')\n",
    "    plt.ylabel('Амплитуда')\n",
    "    plt.grid()\n",
    "    plt.show()\n",
    "\n",
    "plotSignal('S140314_161345')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Данной реализации соответствует 6 проходов, для решения задачи классификации <b>требуется</b>:\n",
    "\n",
    " 1. выделить соответсвующие проходы; \n",
    " 2. для проходов определить (расчитать) признаки;\n",
    " 3. для сформированного признакового пространства добавить вектор меток.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### 1. Выделение соответствующих проходов"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Определим, что длительность прохода должна быть равной 8 секунд (<b>`step_time`</b>), то есть 4000 отсчетов. Визуализируем все реализации и оценим возможность обеспечить длительность прохода 8 с. \n",
    "\n",
    "В ходе рассмотрения установлено, что для корректного автоматизированного выделения проходов (<b>`selection_steps`</b>) требуется дополнительно модифицировать исходные данные, в части точного определения начала прохода (<b>`dop_group_steps`</b>): "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "step_time = 8\n",
    "len_step = step_time * Fs\n",
    "print('Длительность прохода: {} отсчетов'.format(len_step))\n",
    "count_steps = {}\n",
    "for name in name_signal:\n",
    "    if(meta_data[meta_data['nameFile'] == name].shape[0] > 0 and\\\n",
    "       meta_data[meta_data['nameFile'] == name].Steps.values[0] > 0):\n",
    "        count_steps[name] = meta_data[meta_data['nameFile'] == name].Steps.values[0]\n",
    "    else:\n",
    "        name_signal.remove(name)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dop_group_steps = {}\n",
    "dop_group_steps['S140507_180359'] = [8, 28, 48]\n",
    "dop_group_steps['S140507_175631'] = [4, 16, 30]\n",
    "dop_group_steps['S140507_175740'] = [3, 13, 23]\n",
    "dop_group_steps['S140507_175116'] = [3, 15, 27]\n",
    "dop_group_steps['S140507_180656'] = [6, 25, 44]\n",
    "dop_group_steps['S140507_175227'] = [5, 19, 33]\n",
    "dop_group_steps['S140507_181133'] = [6, 18, 28]\n",
    "dop_group_steps['S140507_175006'] = [4, 17, 43]\n",
    "dop_group_steps['S140507_175906'] = [6, 23, 41]\n",
    "dop_group_steps['S140507_175445'] = [4, 16, 30]\n",
    "dop_group_steps['S140507_180527'] = [6, 26, 43]\n",
    "dop_group_steps['S140507_175536'] = [3, 13, 24]\n",
    "dop_group_steps['S140507_180929'] = [4, 15, 26]\n",
    "dop_group_steps['S140507_180021'] = [9, 26, 45]\n",
    "dop_group_steps['S140507_180224'] = [8, 24, 44]\n",
    "dop_group_steps['S140507_181031'] = [4, 17, 27]\n",
    "dop_group_steps['S140507_170801'] = [8, 22, 48, 65]\n",
    "dop_group_steps['S140507_172302'] = [9, 28, 50, 66]\n",
    "dop_group_steps['S140507_171954'] = [3, 12, 25, 37]\n",
    "dop_group_steps['S140507_170942'] = [8, 22, 42, 57]\n",
    "dop_group_steps['S140507_174541'] = [8, 18, 35, 45]\n",
    "dop_group_steps['S140507_181239'] = [4, 17, 30]\n",
    "dop_group_steps['S140507_174732'] = []\n",
    "dop_group_steps['S140507_171831'] = [3, 12, 25, 36]\n",
    "dop_group_steps['S140507_171440'] = []\n",
    "dop_group_steps['S140507_170319'] = [10, 25, 48, 68]\n",
    "dop_group_steps['S140715_163612'] = []\n",
    "dop_group_steps['S140507_178441'] = []\n",
    "dop_group_steps['S140507_165507'] = [8, 28, 52, 72]\n",
    "dop_group_steps['S140507_171708'] = [3, 14, 26, 37]\n",
    "dop_group_steps['S140507_173900'] = [6, 21, 54, 75]\n",
    "dop_group_steps['S140507_175999'] = []\n",
    "dop_group_steps['S140507_173020'] = [2, 13, 24, 35]\n",
    "dop_group_steps['S140507_171330'] = [2, 14, 27, 39]\n",
    "dop_group_steps['S140507_171605'] = [1, 12, 25, 35]\n",
    "dop_group_steps['S140507_170631'] = [6, 25.5, 48.5, 65.5]\n",
    "dop_group_steps['S140507_171100'] = [7, 24, 46, 54.5]\n",
    "dop_group_steps['S140507_172142'] = [8, 20, 44, 62]\n",
    "dop_group_steps['S140507_165330'] = []\n",
    "dop_group_steps['S140507_165315'] = [9, 30, 54, 72]\n",
    "dop_group_steps['S140507_174319'] = [4, 15 ,33, 45]\n",
    "dop_group_steps['S140507_165642'] = [10, 29, 52, 70.5]\n",
    "dop_group_steps['S140507_174435'] = [5, 15, 33, 44]\n",
    "dop_group_steps['S140507_172534'] = [4, 15, 26, 38]\n",
    "dop_group_steps['S140507_172530'] = []\n",
    "dop_group_steps['S140507_172832'] = [2, 13, 23, 35]\n",
    "dop_group_steps['S140507_172731'] = [2, 12, 24, 35]\n",
    "dop_group_steps['S140507_174039'] = [7, 24, 43, 57]\n",
    "dop_group_steps['S140507_170459'] = [10, 26, 53, 68]\n",
    "dop_group_steps['S140507_172926'] = [2, 12, 25, 35]\n",
    "dop_group_steps['S140507_170006'] = [8, 31, 56, 74]\n",
    "dop_group_steps['S140507_174213'] = [5, 15, 35, 46]\n",
    "\n",
    "correct_group_count_step = dop_group_steps.keys()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def selection_steps(j,\n",
    "                     signal,  \n",
    "                     len_steps, \n",
    "                     count_steps, \n",
    "                     len_step):\n",
    "    ind_low_steps = []\n",
    "\n",
    "    le_ = 0\n",
    "    fla = False\n",
    "    for i in range(count_steps):\n",
    "        \n",
    "        max_val = np.max(signal[i*len_steps:(i+1)*len_steps])\n",
    "        max_index = signal.index(max_val)\n",
    "        low_index = int(max_index-len_step/2)\n",
    "        high_index = int(max_index+len_step/2)\n",
    "        if(low_index < 0):\n",
    "            high_index = high_index - low_index\n",
    "            low_index = 0\n",
    "        sig_s = signal[low_index:high_index]\n",
    "        if(sig_s.index(np.max(sig_s)) != len_step):\n",
    "            delta = int(len_step / 2 - sig_s.index(np.max(sig_s)))\n",
    "            low_index = low_index  - delta\n",
    "            high_index = high_index - delta\n",
    "\n",
    "        if(low_index < 0):\n",
    "            high_index = high_index - low_index\n",
    "            low_index = 0\n",
    "        \n",
    "        ind_low_steps.append(int(low_index)/500)\n",
    "    return ind_low_steps"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for name in name_signal:\n",
    "    if(name not in correct_group_count_step):\n",
    "        step = selection_steps(name, seism[name],\n",
    "                                int(len(seism[name])/count_steps[name]), \n",
    "                                count_steps[name], len_step)\n",
    "        dop_group_steps[name] = step\n",
    "        \n",
    "print ('Количество реализаций для которых определено начало прохода - {}.'.format(len(dop_group_steps)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Соберем все данные в один класс и сформируем словарь экземпляров класса `my_signal` характеризующий  реализации сейсмических сигналов.\n",
    "\n",
    "Обработаем сформированные данные: \n",
    " - центрируем сигнал;\n",
    " - применяем фильтр (ФНЧ и ПФ);\n",
    " - нормируем сигнал."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "class my_signal(object):\n",
    "    \"\"\"docstring\"\"\"\n",
    " \n",
    "    def __init__(self, name, ind, ind_step, signal, targetCountN, fs, typeS):\n",
    "        \"\"\"Constructor\"\"\"\n",
    "        self.name = name\n",
    "        self.ind = ind\n",
    "        self.ind_step = ind_step\n",
    "        self.signal = signal\n",
    "        self.flagN = -1\n",
    "        self.step_time = -100500\n",
    "        self.pause_time = -100500\n",
    "        self.step_count = -100500\n",
    "        self.period = -100500\n",
    "        self.energy = -100500\n",
    "        self.countPolin = -100500\n",
    "        self.centr = 1\n",
    "        if(targetCountN != -1):\n",
    "            if(targetCountN == 1):\n",
    "                self.targetCountN = u'Одиночный'\n",
    "                self.flagN = 0\n",
    "            else:\n",
    "                self.targetCountN = u'Групповой'\n",
    "                self.flagN = 1\n",
    "        else:\n",
    "            self.targetCountN = ''\n",
    "            self.flagN = -1\n",
    "        self.fs = fs\n",
    "        self.typeS = typeS\n",
    "    def info(self):\n",
    "        print('--------------------------------')\n",
    "        print('Название сигнала: {}'.format(self.name))\n",
    "        print('Индекс сигнала в выборке: {}'.format(self.ind))\n",
    "        print('Номер прохода в реализации: {}'.format(self.ind_step))\n",
    "        print('Объект которому соответствует сигнал: {}'.format(self.targetCountN))\n",
    "        print('Частота дискретизации сигнала: {}'.format(self.fs))\n",
    "    # Визуализация сигнала во временной области    \n",
    "    def plot_signal(self, \n",
    "                    size = (15, 6), \n",
    "                    delTime = 2,\n",
    "                    nameFile = '',\n",
    "                    dpi2 = 250,\n",
    "                    color = 'r'):\n",
    "        time_x = np.linspace(0,len(self.signal)/self.fs, len(self.signal))\n",
    "        plt.title(u'Наименование :' + self.name + '\\n ' + str(self.targetCountN) + ' ' + '\\n')\n",
    "        plt.plot(time_x,self.signal, color)\n",
    "        plt.legend([self.typeS])\n",
    "        plt.xlabel(u'Время, с')\n",
    "        plt.ylabel(u'Амплитуда')\n",
    "        plt.grid() \n",
    "        plt.savefig(nameFile, dpi = 250)\n",
    "    # Визуализация спекта сигнала        \n",
    "    def plot_spectrum(self, size = (15, 6), delFreq = 50):\n",
    "        fig = plt.figure(figsize = size)\n",
    "        plt.plot(self.frq,abs(self.spectr),'r')\n",
    "        plt.title(u'Спектр сигнала')\n",
    "        plt.xlabel(u'Частота (Гц)')\n",
    "        plt.ylabel(u'|Y(freq)|')\n",
    "        ax = fig.gca()\n",
    "        ax.set_xticks(numpy.arange(0, int(self.fs/2), delFreq))\n",
    "        plt.grid()   \n",
    "    # Вычисление спекта сигнала            \n",
    "    def creat_spectrum(self):\n",
    "        n = len(self.signal) # Длительность сигнала\n",
    "        k = arange(n)\n",
    "        T = n/self.fs\n",
    "        frq = k/T \n",
    "        frq = frq[range(int(n/2))] # Диапазон частот\n",
    "\n",
    "        Y = np.fft(self.signal)/n # Вычисление быстрого преобразования Фурье и его нормализация\n",
    "        Y = np.abs(Y[range(int(n/2))])\n",
    "        self.spectr = Y\n",
    "        self.frq = frq\n",
    "    # Фильтрация сигнала (ФНЧ, ФВЧ, полосовой фильтр и режекторный фильтр)             \n",
    "    def signal_filter(self, cutoff, order=5, btypeFilter = 'low'):\n",
    "        self2 = copy.copy(self)\n",
    "        nyq = 0.5 * self.fs\n",
    "        if((btypeFilter == 'low') or (btypeFilter == 'highpass')):\n",
    "            normal_cutoff = cutoff / nyq\n",
    "            b, a = signal.butter(order, normal_cutoff, btype = btypeFilter, analog=False)\n",
    "        elif((btypeFilter == 'bandpass') or (btypeFilter == 'bandstop')):\n",
    "            normal_cutoff = [cutoff[0]/nyq, cutoff[1]/nyq]\n",
    "            b, a = signal.butter(order, normal_cutoff, btype=btypeFilter, analog=False)\n",
    "        self2.signal = signal.lfilter(b, a, self.signal)\n",
    "        return self2\n",
    "\n",
    "    def operation_erase_mean(self):\n",
    "        me = np.mean(self.signal) \n",
    "        self.signal = self.signal - me\n",
    "    # Нормировка сигнала    \n",
    "    def operation_norm(self):\n",
    "        self.signal = self.signal / np.max(np.abs(self.signal))\n",
    "        \n",
    "    def operation_high_low(self):\n",
    "        \n",
    "\n",
    "        p75 = np.percentile(self.signal, 75)\n",
    "        p25 = np.percentile(self.signal, 25)\n",
    "\n",
    "        qr = p75 - p25\n",
    "\n",
    "        self.high = p75 + 1.5 * qr\n",
    "        self.low = p25 - 1.5 * qr\n",
    "        \n",
    "    def operation_set_high_low(self, high, low):\n",
    "        self.high = high\n",
    "        self.low = low "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "\n",
    "big_seism = {}\n",
    "\n",
    "for name in dop_group_steps:\n",
    "    if(len(dop_group_steps[name]) > 0):\n",
    "        temp_meta = meta_data[meta_data['nameFile'] == name]\n",
    "        sig_s_big = my_signal(name=name,\n",
    "                         ind = temp_meta.index.values[0],\n",
    "                         ind_step = 0,\n",
    "                         signal = seism[name],\n",
    "                         targetCountN = temp_meta.Count.values[0],\n",
    "                         fs = Fs,\n",
    "                         typeS = 'Сейсмический сигнал')\n",
    "        sig_s_big = sig_s_big.signal_filter(btypeFilter='low', cutoff= 35)\n",
    "        sig_s_big = sig_s_big.signal_filter(btypeFilter='bandstop', cutoff= [45, 55])\n",
    "#         sig_s_big.operation_norm()\n",
    "        sig_s_big.operation_high_low()\n",
    "        big_seism[name]=sig_s_big\n",
    "\n",
    "print ('Количество \"больших\" сигналов - {} реализаций.\\n'.format(len(big_seism)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Для примера визуализируем два экземпляра для первого и второго класса:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(figsize = (15, 6))\n",
    "subplot(1,2,1)\n",
    "big_seism['S140314_161345'].plot_signal()\n",
    "subplot(1,2,2)\n",
    "big_seism['S140507_165315'].plot_signal()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Создадим новый словарь экземпляров класса `my_signal` для \"малых\" сигналов - проходов."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "\n",
    "small_seism = []\n",
    "\n",
    "for (j,name) in enumerate(big_seism.keys()):\n",
    "    sig = big_seism[name]\n",
    "\n",
    "    temp_meta = meta_data[meta_data['nameFile'] == sig.name]\n",
    "\n",
    "    for (i, ind) in enumerate(dop_group_steps[sig.name]):\n",
    "\n",
    "        if(len(dop_group_steps[sig.name]) > 0):\n",
    "            sig_s_small = my_signal(name=sig.name,\n",
    "                             ind = sig.ind,\n",
    "                             ind_step = i,\n",
    "                             signal = sig.signal[int(ind*sig.fs):int((ind+8)*sig.fs)],\n",
    "                             targetCountN = temp_meta.Count.values[0],\n",
    "                             fs = sig.fs,\n",
    "                             typeS = 'Сейсмический сигнал')\n",
    "            \n",
    "            sig_s_small.operation_erase_mean()\n",
    "            \n",
    "\n",
    "            sig_s_small.operation_set_high_low(high=sig.high, low=sig.low)\n",
    "            #sig_s_small.operation_high_low()\n",
    "            #sig_s_small.feature_time(low_time = 50)\n",
    "\n",
    "            \n",
    "            small_seism.append(sig_s_small)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print ('Количество \"малых\" сигналов - {} проходов.\\n'.format(len(small_seism)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(figsize = (15, 6))\n",
    "small_seism[800].info()\n",
    "small_seism[800].plot_signal()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "##### 2. Для проходов определим (расчитаем) признаки:\n",
    "\n",
    " 1. средняя длительность фазы переноса - период переноса ноги;\n",
    " 2. средняя длительность фазы опоры - постановка ноги на опору;\n",
    " 3. средняя длительность полного цикла (период двойного шага) - слагается для каждой ноги из фазы опоры и фазы переноса конечности;\n",
    " 4. энергия сигнала;\n",
    " 5. количество шагов.\n",
    "\n",
    "Для расчета 1,2,3 и 5 признака требуется выделить из прохода фазу переноса или фазу опоры.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "indEx = 800\n",
    "seismic_signal = small_seism[indEx].signal\n",
    "N = len(seismic_signal)\n",
    "step_pause = [0] * len(seismic_signal)\n",
    "\n",
    "st = [0] * N\n",
    "\n",
    "p75 = np.percentile(seismic_signal, 75)\n",
    "p25 = np.percentile(seismic_signal, 25)\n",
    "\n",
    "qr = p75 - p25\n",
    "high = p75 + 1.5 * qr\n",
    "low = p25 - 1.5 * qr\n",
    "\n",
    "for j in range(len(seismic_signal)-1):\n",
    "    i = j + 1\n",
    "    if(seismic_signal[i] >= high or seismic_signal[i] <= low):\n",
    "        st[i] = 1\n",
    "        step_pause[i] = seismic_signal[i]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig = plt.figure(figsize = (15, 6))\n",
    "stepsS = pd.DataFrame()\n",
    "for i in range(len(st)-1):\n",
    "    if(st[i] == 0 and st[i+1] == 1 and len(stepsS) > 0):\n",
    "        stepsS.set_value(len(stepsS)-1,'finish', (i))\n",
    "    elif(st[i] == 1 and st[i+1] == 0):\n",
    "        stepsS.set_value(len(stepsS),'start', i)\n",
    "stepsS = stepsS[:-1]\n",
    "stepsS['time_pause'] = stepsS['finish'] - stepsS['start']\n",
    "stepsS['time_pause'].plot()\n",
    "stepsS = stepsS[(stepsS['time_pause']>50)]\n",
    "stepsS['time_pause'].plot()\n",
    "plt.legend(['Фазы переноса до фильтрации', 'Фазы переноса после фильтрации'])\n",
    "plt.xlabel(u'Номер фазы переноса')\n",
    "plt.ylabel(u'Время * 500 Гц')\n",
    "stepsS['time_pause'].median()\n",
    "plt.grid()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "walk = []\n",
    "fig = plt.figure(figsize = (15, 6))\n",
    "#plt.plot(seismic_signal_by_filter_30)\n",
    "time_x = np.linspace(0,len(seismic_signal)/Fs, len(seismic_signal))\n",
    "plt.plot(time_x, seismic_signal)\n",
    "for i in range(stepsS.shape[0]-1):\n",
    "    plt.axvline(stepsS.iloc[i]['finish']/Fs, color = 'red')\n",
    "    plt.axvline(stepsS.iloc[i+1]['start']/Fs, color = 'green')\n",
    "    \n",
    "    stepsS.set_value(stepsS.index[i],'time_walk',stepsS.iloc[i+1]['start'] - stepsS.iloc[i]['finish'])\n",
    "plt.legend(['Расматриваемый сигнал', 'Отсечка начала фазы опоры', 'Отсечка окончания фазы опоры'])\n",
    "plt.xlabel(u'Время, с')\n",
    "plt.ylabel(u'Амплитуда')\n",
    "plt.grid()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Таким образом разработан механизм выделения фазы переноса и опоры. \n",
    "\n",
    "Добавим в сформированный класс метод по выделению признаков -  <b>`feature_time`</b>."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "class my_signal(object):\n",
    "    \"\"\"docstring\"\"\"\n",
    " \n",
    "    def __init__(self, name, ind, ind_step, signal, targetCountN, fs, typeS):\n",
    "        \"\"\"Constructor\"\"\"\n",
    "        self.name = name\n",
    "        self.ind = ind\n",
    "        self.ind_step = ind_step\n",
    "        self.signal = signal\n",
    "        self.flagN = -1\n",
    "        self.step_time = -100500\n",
    "        self.pause_time = -100500\n",
    "        self.step_count = -100500\n",
    "        self.period = -100500\n",
    "        self.energy = -100500\n",
    "        self.countPolin = -100500\n",
    "        self.centr = 1\n",
    "        if(targetCountN != -1):\n",
    "            if(targetCountN == 1):\n",
    "                self.targetCountN = u'Одиночный'\n",
    "                self.flagN = 0\n",
    "            else:\n",
    "                self.targetCountN = u'Групповой'\n",
    "                self.flagN = 1\n",
    "        else:\n",
    "            self.targetCountN = ''\n",
    "            self.flagN = -1\n",
    "        self.fs = fs\n",
    "        self.typeS = typeS\n",
    "    def info(self):\n",
    "        print('--------------------------------')\n",
    "        print('Название сигнала: {}'.format(self.name))\n",
    "        print('Индекс сигнала в выборке: {}'.format(self.ind))\n",
    "        print('Номер прохода в реализации: {}'.format(self.ind_step))\n",
    "        print('Объект которому соответствует сигнал: {}'.format(self.targetCountN))\n",
    "        print('Частота дискретизации сигнала: {}'.format(self.fs))\n",
    "        \n",
    "    def plot_signal(self, \n",
    "                    size = (15, 6), \n",
    "                    delTime = 2,\n",
    "                    nameFile = '',\n",
    "                    dpi2 = 250,\n",
    "                    color = 'r'):\n",
    "        time_x = np.linspace(0,len(self.signal)/self.fs, len(self.signal))\n",
    "        plt.title(u'Наименование :' + self.name + '\\n ' + str(self.targetCountN) + ' ' + '\\n')\n",
    "        plt.plot(time_x,self.signal, color)\n",
    "        plt.legend([self.typeS])\n",
    "        plt.xlabel(u'Время, с')\n",
    "        plt.ylabel(u'Амплитуда')\n",
    "\n",
    "        plt.savefig(nameFile, dpi = 250)\n",
    "            \n",
    "    def plot_spectrum(self, size = (15, 6), delFreq = 50):\n",
    "        fig = plt.figure(figsize = size)\n",
    "        plt.plot(self.frq,abs(self.spectr),'r')\n",
    "        plt.title(u'Спектр сигнала')\n",
    "        plt.xlabel(u'Частота (Гц)')\n",
    "        plt.ylabel(u'|Y(freq)|')\n",
    "        ax = fig.gca()\n",
    "        ax.set_xticks(numpy.arange(0, int(self.fs/2), delFreq))\n",
    "        plt.grid()   \n",
    "        \n",
    "    def creat_spectrum(self):\n",
    "        n = len(self.signal) # Длительность сигнала\n",
    "        k = arange(n)\n",
    "        T = n/self.fs\n",
    "        frq = k/T \n",
    "        frq = frq[range(int(n/2))] # Диапазон частот\n",
    "\n",
    "        Y = np.fft(self.signal)/n # Вычисление быстрого преобразования Фурье и его нормализация\n",
    "        Y = np.abs(Y[range(int(n/2))])\n",
    "        self.spectr = Y\n",
    "        self.frq = frq\n",
    "        \n",
    "    def signal_filter(self, cutoff, order=5, btypeFilter = 'low'):\n",
    "        self2 = copy.copy(self)\n",
    "        nyq = 0.5 * self.fs\n",
    "        if((btypeFilter == 'low') or (btypeFilter == 'highpass')):\n",
    "            normal_cutoff = cutoff / nyq\n",
    "            b, a = signal.butter(order, normal_cutoff, btype = btypeFilter, analog=False)\n",
    "        elif((btypeFilter == 'bandpass') or (btypeFilter == 'bandstop')):\n",
    "            normal_cutoff = [cutoff[0]/nyq, cutoff[1]/nyq]\n",
    "            b, a = signal.butter(order, normal_cutoff, btype=btypeFilter, analog=False)\n",
    "        self2.signal = signal.lfilter(b, a, self.signal)\n",
    "        return self2\n",
    "    \n",
    "    def operation_erase_mean(self):\n",
    "        me = np.mean(self.signal) \n",
    "        self.signal = self.signal - me\n",
    "        \n",
    "    def operation_norm(self):\n",
    "        self.signal = self.signal / np.max(np.abs(self.signal))\n",
    "        \n",
    "    def print_features(self):\n",
    "        \n",
    "        print ('Средняя длительность шага:{}'.format(self.step_time))\n",
    "        print ('Средняя длительность паузы:{}'.format(self.pause_time))\n",
    "        print ('Количество шагов:{}'.format(self.step_count))\n",
    "        print ('Энергия сигнала:{}'.format(self.energy))\n",
    "    def operation_high_low(self):\n",
    "        \n",
    "\n",
    "        p75 = np.percentile(self.signal, 75)\n",
    "        p25 = np.percentile(self.signal, 25)\n",
    "\n",
    "        qr = p75 - p25\n",
    "\n",
    "        self.high = p75 + 1.5 * qr\n",
    "        self.low = p25 - 1.5 * qr\n",
    "        \n",
    "    def operation_set_high_low(self, high, low):\n",
    "        self.high = high\n",
    "        self.low = low \n",
    "    def feature_time(self, low_time):\n",
    "        if(self.flagN != -1):\n",
    "            \n",
    "            N = len(self.signal)\n",
    "\n",
    "            signal2 = self.signal\n",
    "\n",
    "            \n",
    "            st = [0] * N\n",
    "\n",
    "            for j in range(N-1):\n",
    "                i = j + 1\n",
    "                if(signal2[i] >= self.high or signal2[i] <= self.low):\n",
    "                    st[i] = 1\n",
    "                else:\n",
    "                    signal2[i] = 0    \n",
    "\n",
    "            start = []\n",
    "            finish = []\n",
    "            for i in range(len(st)-1):\n",
    "                if(st[i] == 0 and st[i+1] == 1 and len(start) > 0):\n",
    "                    finish.append(i)\n",
    "                elif(st[i] == 1 and st[i+1] == 0):\n",
    "                    start.append(i)\n",
    "\n",
    "            pause_time_2 = []\n",
    "            finish_2 = []\n",
    "            start_2 = []\n",
    "            step_time_2 = []\n",
    "            pause_time = list(map(lambda x: x[0] - x[1], zip(finish, start)))\n",
    "\n",
    "            for i in range(len(pause_time)):\n",
    "                if(pause_time[i] > low_time):\n",
    "                    pause_time_2.append(pause_time[i])\n",
    "                    finish_2.append(finish[i])\n",
    "                    start_2.append(start[i])\n",
    "\n",
    "            for i in range(len(pause_time_2)-1):\n",
    "                step_time_2.append(start_2[i+1] - finish_2[i])\n",
    "\n",
    "            pause_time = []\n",
    "            finish = []\n",
    "            start = []\n",
    "            step_time = []\n",
    "            period = []\n",
    "            for i in range(len(pause_time_2)-1):\n",
    "                if(pause_time_2[i] < 5000 and step_time_2[i] < 400 and step_time_2[i] > 2 ):\n",
    "                    pause_time.append(pause_time_2[i])\n",
    "                    finish.append(finish_2[i])\n",
    "                    start.append(start_2[i])   \n",
    "                    step_time.append(step_time_2[i])\n",
    "                    period.append(step_time_2[i] + pause_time_2[i])\n",
    "\n",
    "            self.step_time = np.mean(step_time)\n",
    "            self.pause_time = np.mean(pause_time)\n",
    "            self.step_count = len(step_time)\n",
    "            self.period = np.mean(period)\n",
    "\n",
    "            analytic_signal = hilbert(signal2)\n",
    "            amplitude_envelope = np.abs(analytic_signal)\n",
    "            signal2 = np.abs(signal2) / np.max(signal2)\n",
    "            self.energy = np.sum(signal2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Для всего словаря сигналов применим метод `feature_time`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%%time\n",
    "\n",
    "small_seism_feature = []\n",
    "for signal in small_seism:\n",
    "    \n",
    "    sig_s_small_feature = my_signal(name=signal.name,\n",
    "                     ind = signal.ind,\n",
    "                     ind_step = signal.ind_step,\n",
    "                     signal = signal.signal,\n",
    "                     targetCountN = signal.targetCountN,\n",
    "                     fs = signal.fs,\n",
    "                     typeS = signal.typeS)\n",
    "\n",
    "    sig_s_small_feature.flagN = signal.flagN\n",
    "    sig_s_small_feature.operation_set_high_low(high=signal.high, low=signal.low)\n",
    "    \n",
    "    sig_s_small_feature.feature_time(low_time = 50)\n",
    "\n",
    "    small_seism_feature.append(sig_s_small_feature)\n",
    "print ('Количество \"малых\" сигналов - {} проходов с извлеченными признаками.\\n'.format(len(small_seism_feature)))    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "##### 3. для сформированного признакового пространства добавить вектор меток.\n",
    "\n",
    "Выделим из словаря сигналов необходимые признаки и сформируем DataFrame"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X = pd.DataFrame()\n",
    "X['step_time'] = list(map(lambda x: x.step_time, small_seism_feature))\n",
    "X['pause_time'] = list(map(lambda x: x.pause_time, small_seism_feature))\n",
    "X['step_count'] = list(map(lambda x: x.step_count, small_seism_feature))\n",
    "X['energy'] = list(map(lambda x: x.energy, small_seism_feature))\n",
    "X['target'] = list(map(lambda x: x.flagN, small_seism_feature))\n",
    "X['period'] = list(map(lambda x: x.period, small_seism_feature))\n",
    "X = X[X['target'] > -1]\n",
    "X = X.dropna()\n",
    "print ('Количество групповых целей в выборке {}'.format(X[X['target'] == 1].shape[0]))\n",
    "print ('Количество одиночных целей в выборке {}'.format(X[X['target'] == 0].shape[0]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X.head()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Визуализируем признаки"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(figsize(20, 5))\n",
    "plt.subplot(1, 5, 1)\n",
    "df = pd.melt(X, value_vars=['pause_time'], id_vars='target')\n",
    "sns.violinplot(x='variable', y ='value', hue='target', data=df, scale ='count', split=True,palette=\"Set1\")\n",
    "plt.grid()\n",
    "plt.subplot(1, 5, 2)\n",
    "df = pd.melt(X, value_vars=['step_time'], id_vars='target')\n",
    "sns.violinplot(x='variable', y ='value', hue='target', data=df, scale ='count', split=True,palette=\"Set1\")\n",
    "plt.grid()\n",
    "plt.subplot(1, 5, 3)\n",
    "df = pd.melt(X, value_vars=['step_count'], id_vars='target')\n",
    "sns.violinplot(x='variable', y ='value', hue='target', data=df, scale ='count', split=True,palette=\"Set1\")\n",
    "plt.grid()\n",
    "plt.subplot(1, 5, 4)\n",
    "df = pd.melt(X, value_vars=['energy'], id_vars='target')\n",
    "sns.violinplot(x='variable', y ='value', hue='target', data=df, scale ='count', split=True,palette=\"Set1\")\n",
    "plt.grid()\n",
    "plt.subplot(1, 5, 5)\n",
    "df = pd.melt(X, value_vars=['period'], id_vars='target')\n",
    "sns.violinplot(x='variable', y ='value', hue='target', data=df, scale ='count', split=True,palette=\"Set1\")\n",
    "plt.grid()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Рассмотрим попарную корреляцию между признаками"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(figsize(6, 5))\n",
    "sns.heatmap(X.corr(), cmap='PuOr');"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Заметим, что:\n",
    " * признаки step_time, pause_time сильно положительно коррелируют с признаком period;\n",
    " * признак pause_time отрицательно коррелирует с признаком energy.\n",
    " \n",
    "Такие коэффициенты корреляции обусловлены природой формирования признакового пространства во временной области.\n",
    "\n",
    "Кроме того, в рамках <b>исследования физиологии движений и физиологии активности</b>, Николай Александрович Бернштейн оценил, что при ходьбе в среднем темпе фаза опоры длится примерно 60 % от цикла двойного шага, фаза переноса примерно 40 %, оценим получившиеся данные. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Средняя длительность двойного шага - {}'.format(np.round(X.period.mean(), 3)))\n",
    "print('Средняя длительность фазы переноса - {} ({} % от цикла двойного шага)'.format(np.round(X.step_time.mean(),3), \n",
    "                                                                                     np.round(X.step_time.mean() / X.period.mean(), 3)))\n",
    "print('Средняя длительность фазы опоры - {} ({} % от цикла двойного шага)'.format(np.round(X.pause_time.mean(),3), \n",
    "                                                                                     np.round(X.pause_time.mean() / X.period.mean(), 3)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Представленные расчеты в целом совпадают с теоретическими результатами"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "y = X['target']\n",
    "X = X.drop('target', 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "forest = RandomForestClassifier(n_estimators=40)\n",
    "\n",
    "forest.fit(X,y)\n",
    "\n",
    "importances = forest.feature_importances_\n",
    "std = np.std([tree.feature_importances_ for tree in forest.estimators_],\n",
    "             axis=0)\n",
    "ind = np.argsort(importances)[::-1]\n",
    "\n",
    "# Print the feature ranking\n",
    "print(\"Сортировка признаков по информативности:\")\n",
    "\n",
    "for f in range(X.shape[1]):\n",
    "    print(\"%d. Признак -  %s (%f)\" % (f + 1, X.columns[ind[f]], importances[ind[f]]))\n",
    "\n",
    "# Plot the feature importances of the forest\n",
    "plt.figure()\n",
    "plt.title(\"Ранг признаков по информативности\")\n",
    "plt.bar(range(X.shape[1]), importances[indices],\n",
    "       yerr=std[indices], align=\"center\")\n",
    "plt.xticks(range(X.shape[1]), indices)\n",
    "plt.grid()\n",
    "plt.xlim([-1, X.shape[1]])\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Кросс-валидация"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Разделим выборку на обучающую и тестовую"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.cross_validation import train_test_split\n",
    "X_train, X_test, y_train, y_test = train_test_split(X, y, \n",
    "                                                    test_size=0.3, \n",
    "                                                    random_state=42,\n",
    "                                                   stratify = y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "forest = RandomForestClassifier(n_estimators=500)\n",
    "forest.fit(X_train, y_train)\n",
    "test_pred = forest.predict(X_test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "accuracy_score(y_test, test_pred)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "y_test.value_counts()[0] / y_test.shape[0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f1_score(y_test, test_pred)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Подберем параметры для модели"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rf_tree_params = {'n_estimators': (50, 100, 150),\n",
    "                  'min_samples_split' : list(range(1,5)),\n",
    "               'max_depth': list(range(6,14))}\n",
    "n_folds = 5"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "grid = GridSearchCV(clf, rf_tree_params, cv=n_folds, \n",
    "                    n_jobs=-1)\n",
    "grid.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "clf_best_score = grid.best_score_\n",
    "clf_best_params = grid.best_params_\n",
    "clf_best = grid.best_estimator_\n",
    "mean_validation_scores = []\n",
    "print(\"Лучший результат\", clf_best_score)\n",
    "print(\"лучшие параметры\", clf_best_params)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Проверяем сходимость модели"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_with_std(x, data, **kwargs):\n",
    "        mu, std = data.mean(1), data.std(1)\n",
    "        lines = plt.plot(x, mu, '-', **kwargs)\n",
    "        plt.fill_between(x, mu - std, mu + std, edgecolor='none',\n",
    "                         facecolor=lines[0].get_color(), alpha=0.2)\n",
    "        \n",
    "def plot_learning_curve(clf, X, y, scoring, cv=5):\n",
    " \n",
    "    train_sizes = np.linspace(0.05, 1, 20)\n",
    "    n_train, val_train, val_test = learning_curve(clf,\n",
    "                                                  X, y, train_sizes, cv=cv,\n",
    "                                                  scoring=scoring)\n",
    "    plot_with_std(n_train, val_train, label='training scores', c='green')\n",
    "    plot_with_std(n_train, val_test, label='validation scores', c='red')\n",
    "    plt.xlabel('Training Set Size'); plt.ylabel(scoring)\n",
    "    plt.legend()\n",
    "    plt.grid()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_learning_curve(RandomForestClassifier(n_estimators=clf_best_params['n_estimators'], \n",
    "                                           max_depth=clf_best_params['max_depth'], \n",
    "                                           min_samples_leaf=clf_best_params['min_samples_leaf']),\n",
    "                   X_train, y_train, scoring='f1', cv=10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plot_learning_curve(RandomForestClassifier(n_estimators=clf_best_params['n_estimators'], \n",
    "                                           max_depth=clf_best_params['max_depth'], \n",
    "                                           min_samples_leaf=clf_best_params['min_samples_leaf']),\n",
    "                   X_train, y_train, scoring='accuracy', cv=10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Рассмотрим как влиет количетсво деревьев в модели "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_validation_curve(clf, X, y, cv_param_name, \n",
    "                          cv_param_values, scoring):\n",
    "\n",
    "    val_train, val_test = validation_curve(clf, X, y, cv_param_name,\n",
    "                                           cv_param_values, cv=5,\n",
    "                                                  scoring=scoring)\n",
    "    plot_with_std(cv_param_values, val_train, \n",
    "                  label='training scores', c='green')\n",
    "    plot_with_std(cv_param_values, val_test, \n",
    "                  label='validation scores', c='red')\n",
    "    plt.xlabel(cv_param_name); plt.ylabel(scoring)\n",
    "    plt.legend()\n",
    "    plt.grid()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "estimators = np.arange(25, 350, 25)\n",
    "plot_validation_curve(RandomForestClassifier(min_samples_leaf=clf_best_params['min_samples_leaf'], \n",
    "                                             max_depth=clf_best_params['max_depth']), X_train, y_train, \n",
    "                      cv_param_name='n_estimators', \n",
    "                      cv_param_values= estimators,\n",
    "                   scoring='f1')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Так же рассмотрим как влияет глубина деревьев на качество модели"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "depth = np.arange(3, 25)\n",
    "plot_validation_curve(RandomForestClassifier(n_estimators=clf_best_params['n_estimators'],\n",
    "                                             min_samples_leaf=clf_best_params['min_samples_leaf']), X_train, y_train, \n",
    "                      cv_param_name='max_depth', \n",
    "                      cv_param_values= depth,\n",
    "                   scoring='f1')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Заметим что:\n",
    " * на качество модели практически не влиется количество используемых деревьев;\n",
    " * нет смысла применять глубину больше 15."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "release_forest = RandomForestClassifier(n_estimators=clf_best_params['n_estimators'], \n",
    "                                           max_depth=clf_best_params['max_depth'], \n",
    "                                           min_samples_leaf=clf_best_params['min_samples_leaf'])\n",
    "release_forest.fit(X_train, y_train)\n",
    "release_pred = release_forest.predict(X_test)\n",
    "print('accuracy = {}'.format(accuracy_score(y_test, release_pred)))\n",
    "print('f1 = {}'.format(f1_score(y_test, release_pred)))\n",
    "print('roc_auc = {}'.format(roc_auc_score(y_test, release_pred)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Выводы "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Рассмотрен полный цикл обработки и генерирования признаков из сейсмических сигналов.\n",
    "\n",
    "Объем выборки составил более 1000 сигналов. Уникальных объектов (в данном случае людей по которым регистрировались сигналы) было более 10.\n",
    "\n",
    "Признаки были сформированы во временной области, при дальнейшем исследовании планируется рассмотреть признаки из частотной области, отдельные гармоники или огибающие спектра.\n",
    "\n",
    "Конечно для более точной классификации требуется больше самих сигналов, а также больше объектов."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}
